{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "% CGRtools Tutorial\n",
    "% Dr. Ramil Nugmanov; Dr. Timur Madzhidov; Ravil Mukhametgaleev\n",
    "% Mar 25, 2019"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Data types and operations with them\n",
    "\n",
    "(c) 2019, Dr. Ramil Nugmanov; Dr. Timur Madzhidov; Ravil Mukhametgaleev\n",
    "\n",
    "Installation instructions of CGRtools package information and tutorial's files see on `https://github.com/cimm-kzn/CGRtools`\n",
    "\n",
    "NOTE: Tutorial should be performed sequentially from the start. Random cell running will lead to unexpected results. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pkg_resources\n",
    "if pkg_resources.get_distribution('CGRtools').version.split('.')[:2] != ['3', '1']:\n",
    "    print('WARNING. Tutorial was tested on 3.1 version of CGRtools')\n",
    "else:\n",
    "    print('Welcome!')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load data for tutorial\n",
    "from pickle import load\n",
    "from traceback import format_exc\n",
    "\n",
    "with open('molecules.dat', 'rb') as f:\n",
    "    molecules = load(f) # list of MoleculeContainer objects\n",
    "with open('reactions.dat', 'rb') as f:\n",
    "    reactions = load(f) # list of ReactionContainer objects\n",
    "\n",
    "m1, m2, m3, m4 = molecules # molecule\n",
    "m7 = m3.copy()\n",
    "m11 = m3.copy()\n",
    "m11.standardize()\n",
    "m7.standardize()\n",
    "r1 = reactions[0] # reaction\n",
    "m5 = r1.reactants[0]\n",
    "m8 = m7.substructure([4, 5, 6, 7, 8, 9], as_view=False)\n",
    "m10 =  r1.products[0].copy()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CGRtools has subpackage *containers* with data structures classes:\n",
    "\n",
    "* *MoleculeContainer* - for molecular structure\n",
    "* *ReactionContainer* - for chemical reaction \n",
    "* *CGRContainer* - for Condensed Graph of Reaction\n",
    "* *QueryContainer* - queries for substructure search in molecules\n",
    "* *QueryCGRContainer* - queries for substructure search in CGRs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CGRtools.containers import * # import all containers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.1. MoleculeContainer\n",
    "Molecules are represented as undirected graphs. Molecules contain *Atom* objects and *Bond* objects.\n",
    "\n",
    "*Atom* objects are represented as dictionary with unique number for each atom as key.  \n",
    "\n",
    "*Bond* objects are stored as sparse matrix with adjacent atoms pair as keys for rows and columns.\n",
    "\n",
    "Hereafter, *atom number* is unique integer used to enumerate atoms in molecule. Please, don't confuse it with element number in Periodic Table, hereafter called *element number*.\n",
    "\n",
    "Methods for molecule handling and the arguments of *MoleculeContainer* are described below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.meta # dictionary for molecule properties storage. For example, DTYPE/DATUM fields of SDF file are read into this dictionary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1 # MoleculeContainer supports depiction and graphic representation in Jupyter notebooks."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.depict() # depiction returns SVG image in format string"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('molecule.svg', 'w') as f: # saving image to SVG file\n",
    "    f.write(m1.depict())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m_copy = m1.copy() # isolated copy of molecule\n",
    "m_copy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(m1) # get number of atoms in molecule\n",
    "# or \n",
    "m1.atoms_count"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.bonds_count # number of bonds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.atoms_numbers # list of atoms numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this method calculates additional atoms attributes: number of neighbors and hybridization. See below for usage\n",
    "m1.reset_query_marks() # by default this attributes are None (for speed-up)\n",
    "m3.reset_query_marks()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following notations are used for hybridization of atoms. Values are given as numbers below (in parenthesis symbols that are used in SMILES-like signatures are shown):\n",
    "\n",
    "* 1 (s) - all bonds of atom are single, i.e. sp3 hybridization\n",
    "* 2 (d) - atom has one double bond and others are single, i.e. sp2 hybridization\n",
    "* 3 (t) - atom has one triple or two double bonds and other are single, i.e. sp hybridization\n",
    "* 4 (a) - atom is in aromatic ring\n",
    "\n",
    "Neighbors and hybridizations atom attributes are required for **substructure operations** and structure standardization. See below"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# iterate over atoms using its numbers\n",
    "list(m1.atoms())  # works the same as dict.items()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# iterate over bonds using adjacent atoms numbers\n",
    "list(m1.bonds())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access to atom by number\n",
    "m1.atom(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    m1.atom(10) # raise error for absent atom numbers\n",
    "except KeyError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# access to bond using adjacent atoms numbers\n",
    "m1.bond(1, 4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    m1.bond(1, 3) # raise error for absent bond\n",
    "except KeyError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Atom objects are dictinary-like classes which store information about:\n",
    "* element\n",
    "* isotope\n",
    "* charge\n",
    "* multiplicity\n",
    "* xyz coordinates\n",
    "\n",
    "Also atoms has methods for data integrity checks and include some internally used data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = m1.atom(1)\n",
    "\n",
    "# access to information\n",
    "a.element # element symbol\n",
    "# or\n",
    "a['element']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.charge # formal charge\n",
    "# or\n",
    "a['charge']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.multiplicity # atom multiplicity. None if not set\n",
    "# or\n",
    "a['multiplicity']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.isotope # atom isotope. Default isotope if not set. Default isotopes are the same as used in InChI notation\n",
    "# or\n",
    "a['isotope']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.x # coordinates\n",
    "a.y\n",
    "a.z\n",
    "# or\n",
    "a['x']\n",
    "a['y']\n",
    "a['z']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.neighbors # Number of neighboring atoms is calculated by reset_query_marks method as shown above. It is read-only."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.hybridization # Atoms hybridization is calculated by reset_query_marks method as shown above. It is read-only."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    a.hybridization = 2 # Not subsettable. Read-only! Thus error is raised.\n",
    "except AttributeError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Atomic attributes are subsettable.\n",
    "\n",
    "CGRtools has integrity checks for verification of changes induced by user"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.isotope = 16\n",
    "# or \n",
    "a['isotope'] = 16"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.flush_cache() # due to caching used for speed-up one needs to reset cache of molecule to observe changes in atoms and bonds \n",
    "m1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "try:\n",
    "    a.isotope = 0 # raise error. Isotope with mass 0 could not exist\n",
    "except ValueError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# bond objects also are dictionary-like classes which store information about bond order\n",
    "b = m1.bond(3, 4)\n",
    "\n",
    "b.order\n",
    "#\n",
    "b['order']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b.order = 1 # order change also possible\n",
    "# or\n",
    "b['order'] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.flush_cache() #after flushing cashe one could see molecule with changes\n",
    "m1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "try:\n",
    "    b['order'] = 0 # error raises. Bond order 0 is invalid. To delete bond use method described below\n",
    "except ValueError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One should to use `delete_bond` method to break bond"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.delete_bond(3, 4)\n",
    "m1                    #Note that cashe flushing is not required here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Method `delete_atom` removes atom from the molecule"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.delete_atom(3)\n",
    "m1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m_copy # copy unchanged!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Atoms* and *bonds* objects can be converted into integer representation that could be used to classify their types.\n",
    "\n",
    "*Atom* type is represented by 21 bit code rounded to 32 bit integer number:\n",
    "\n",
    "* 7 bits stand for atom number (2 ** 7 - 1 == 127, currently 118 elements are presented in Periodic Table)\n",
    "* 9 bits are used for isotope (511 posibilities, highest known isotope is ~300)\n",
    "* 3 bits stand for formal charge. Charges range from -3 to +3 rescaled to range 0-6\n",
    "* 2 bits are used for multiplicity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "int(a)\n",
    "# 131596 == 0001000 000010000 011 00\n",
    "# 0001000 == 8 Oxygen\n",
    "# 000010000 == 16 isotope\n",
    "# 011 == 3 (3 - 3 = 0) uncharged\n",
    "# 00 == 0 hasn't multiplicity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "int(b)  # bonds are encoded by their order"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(m1.atom_implicit_h(1)) # get number of implicit hydrogens on atom 1\n",
    "print(m1.atom_explicit_h(1)) # get number of explicit hydrogens on atom 1\n",
    "print(m1.atom_total_h(1)) # get total number of hydrogens on atom 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m1.check_valence() # return list of numbers of atoms with invalid valences"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m4 # molecule with valence errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m4.check_valence()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m3.sssr # Method for application of Smallest Set of Smallest Rings algorithm for rings \n",
    "        # identification. Returns list of lists of atoms forming smallest rings"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Connected components.\n",
    "Sometimes molecules has disconnected components (salts etc).\n",
    "\n",
    "One can find them and split molecule to separate components."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "m2 # it's a salt represented as one graph"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m2.connected_components # list of lists of atoms belonging to graph components"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "anion, cation = m2.split() # split molecule to components"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "anion # graph of only one salt component"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cation # graph of only one salt component"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Union of molecules\n",
    "Sometimes it is more convenient to represent salts as ion pair. Otherwise ambiguity could be introduced, for example in reaction of salt exchange:\n",
    "\n",
    "**Ag+** + **NO3-** + **Na+** + **Br-** = **Ag+** + **Br-** + **Na+** + **NO3-**. Reactants and products sets are the same. \n",
    "\n",
    "In this case one can combine anion-cation pair into single graph. It could be convenient way to represent other molecule mixtures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "salt = anion | cation \n",
    "\n",
    "# or \n",
    "\n",
    "anion.union(cation)\n",
    "\n",
    "salt # this graph has disconnected components, it is considered single compound now"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Substructures could be extracted from molecules.\n",
    "\n",
    "By default, returned substructures are read-only projections of original molecule (except attributes of atoms/bonds). \n",
    "\n",
    "Changes in original molecule (bond breaking/formation, atom insertion/deletion, atom/bond attributes changes) will be mirrored in projection.  \n",
    "\n",
    "Projections share the same neighbors and hybridization attributes as in initial molecule even if could be wrong for substructure "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proj = m3.substructure([4,5,6,7,8,9])  # projection of substructure\n",
    "proj"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m3.atom(4).neighbors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proj.atom(4).neighbors # same as in original molecule and not as it should be in substructure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from networkx.exception import NetworkXError\n",
    "try:\n",
    "    proj.reset_query_marks() # change of structure for projections is blocked\n",
    "except NetworkXError:\n",
    "    print(format_exc())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "benzene = m3.substructure([4,5,6,7,8,9], as_view=False) # Substructure could be extracted as isolated graph (not projection)\n",
    "benzene"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "benzene.atom(4).neighbors is None # empty attribute. Substructure is a new molecule here. We need to call reset_query_marks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "benzene.reset_query_marks()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "benzene.atom(4).neighbors # now number of neighbors for atom 4 is 2. It is not 3 as above where projection was used."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note:\n",
    "\n",
    "* Projection of projection also projection of original molecule\n",
    "* Projection can be converted to isolated molecule by calling method copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proj_copy = proj.copy() # turning projection into molecule using \"copy\" method\n",
    "proj_copy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proj_copy.reset_query_marks() # This not a projection any more but a new molecule"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Changes in projection are mirrored. See example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m3.delete_bond(4, 5) # we've deleted bond\n",
    "m3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proj.flush_cache() # remove cached image. Note that in projection the bond is also deleted.\n",
    "proj"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`augmented_substructure` is a substructure consisting from atoms and a given number of shells of neighboring atoms around it.\n",
    "**deep** argument is a number of considered shells. \n",
    "\n",
    "It also returns projection by default."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "aug = m3.augmented_substructure([10], deep=2) #  atom 10 is Nitrogen\n",
    "aug"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "aug.atom(10).hybridization #atom has two incident double bond."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Atoms Ordering.\n",
    "This functionality is used for canonic numbering of atoms in molecules. Prime number multiplication based Morgan algorithm is used for atom ranking. Property `atoms_order` returns dictionary of atom numbers as keys and their ranks according to canonicalization as values. Equal rank mean that atoms are symmetric (are mapped to each other in automorhisms). In present version, instead of sequential ranks prime numbers are returned.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m5.atoms_order"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Atom number can be changed by `remap` method.\n",
    "\n",
    "This method is useful when it is needed to change order of atoms in molecules. First argument to `remap` method is dictionary with existing atom numbers as keys and desired atom number as values. It is possible to change atom numbers for only part of atoms. Atom numbers could be non-sequencial but need to be unique. \n",
    "\n",
    "If argument *copy* is set ***True*** new object will be created, else existing molecule will be changed. Default is ***False***. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for n, a in m5.atoms():\n",
    "    print(n, a.element)\n",
    "for n, m, b in m5.bonds():\n",
    "    print(m, n, b.order)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "remapped = m5.remap({4:2}, copy=True)\n",
    "remapped"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for n, a in remapped.atoms():\n",
    "    print(n, a.element)\n",
    "for n, m, b in remapped.bonds():\n",
    "    print(m, n, b.order)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.2. ReactionContainer\n",
    "\n",
    "*ReactionContainer* objects has the following properties:\n",
    "\n",
    "* **reactants** - list of reactants molecules\n",
    "* **reagents** - list of reagents molecules\n",
    "* **products** - list of products molecules\n",
    "* **meta** - dictinary of reaction metadata (DTYPE/DATUM block in RDF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r1 # depiction supported"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "r1.meta"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(r1.reactants, r1.products)  # Access to lists of reactant and products. Molecules' signatures are returned by print() method.\n",
    "reactant1, reactant2, reactant3 = r1.reactants\n",
    "product = r1.products[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Reactions also has `standardize`, `aromatize`, `reset_query_marks`, `implicify_hydrogens` and `explicify_hydrogens` methods (see part 3). These methods are applied independently to every molecule in reaction."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.3. CGR\n",
    "*CGRContainer* object is similar to *MoleculeConrtainer*, except some methods. The following methods are not suppoted for *CGRContainer*:\n",
    "\n",
    "* aromatize\n",
    "* standardize\n",
    "* implicify_hydrogens\n",
    "* explicify_hydrogens\n",
    "* atom_implicit_h\n",
    "* atom_explicit_h\n",
    "* atom_total_h\n",
    "* check_valence\n",
    "\n",
    "*CGRContainer* also has some methods absent in *MoleculeConrtainer*:\n",
    "\n",
    "* centers_list\n",
    "* center_atoms\n",
    "* center_bonds\n",
    "\n",
    "*CGRContainer* is undirected graph. Atoms and bonds in CGR has two states: reactant and product."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Composing to CGR\n",
    "\n",
    "As mentioned above, atoms in *MoleculeContainer* have unique numbers. These numbers are used as atom-to-atom mapping in CGRtools upon CGR creation. Thus, atom order for molecules in reaction should correspond to atom-to-atom mapping.  \n",
    "\n",
    "Pair of molecules can be transformed into CGR. Notice that, the same atom numbers in reagents and products imply the same atoms.\n",
    "\n",
    "Reaction also can be composed into CGR. Atom numbers of molecules in reaction are used as atom-to-atom mapping of reactants to products."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cgr1 = m7 ^ m8 # CGR from molecules\n",
    "# or \n",
    "m7.compose(m8)\n",
    "\n",
    "print(cgr1)  # CGRcontainer object currently doesn't support depiction. CGR signature is printed out."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is CGR (depiction is made by ChemAxon externally). You can see changed bonds connected to ring.\n",
    "![cgr1.png](cgr1.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cgr2 = ~r1 # CGR from reactions\n",
    "\n",
    "# or \n",
    "\n",
    "r1.compose()\n",
    "print(cgr2) # signature is printed out."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is CGR for reaction (depiction is made by ChemAxon externally).\n",
    "![cgr2.png](cgr2.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cgr2.reset_query_marks() # CGRs also has reset_query_marks method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = cgr2.atom(2) # atom access is the same as for MoleculeContainer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.element # element attribute\n",
    "\n",
    "# or \n",
    "\n",
    "a['element']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.isotope # isotope attribute\n",
    "\n",
    "# or \n",
    "\n",
    "a['isotope']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For *CGRContainer* attributes `charge`, `multiplicity`, `x`, `y`, `z`, `neighbors` and `hybridization` refer to atom state in reactant of reaction;  arguments `p_charge`, `p_multiplicity`, `p_x`, `p_y`, `p_z`, `p_neighbors` and `p_hybridization` could be used to extract atom state in product part in reaction."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.charge # charge of atom in reactant\n",
    "\n",
    "# or \n",
    "\n",
    "a['charge']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.p_charge # charge of atom in product\n",
    "\n",
    "#or \n",
    "\n",
    "a['p_charge']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.p_multiplicity # multiplicity of atom in product. It is None and thus not returned"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.p_x  # coordinates of atom in product\n",
    "a.p_y\n",
    "a.p_z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.neighbors # number of neighbors of atom in reactant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.p_neighbors # number of neighbors of atom in product"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.hybridization # hybridization of atom in reactant. 1 means only single bonds are incident to atom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a.p_hybridization # hybridization of atom in product. 1 means only single bonds are incident to atom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = cgr1.bond(4, 10) # take bond"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Bonds has `order` and `p_order` attribute\n",
    "If `order` attribute value is ***None***, it means that bond was formed  \n",
    "If `p_order` is ***None***, it means that bond was broken  \n",
    "\n",
    "Both `order` and `p_order` can't be ***None***"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b.order # bond order in reactant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b.p_order is None # bond order in product in None"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**CGR can be decomposed back to reaction**, i.e. reactants and products.\n",
    "\n",
    "Notice that CGR can lose information in case of unbalanced reactions (where some atoms of reactant does not have counterpart in product, and vice versa). Decomposition of CGRs for unbalanced reactions back to reaction may lead to strange (and erroneous) structures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reactant_part, product_part = ~cgr1 # CGR of unbalanced reaction is decomposed back into reaction\n",
    "# or \n",
    "cgr1.decompose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reactant_part # reactants extracted. One can notice it is initial molecule"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "product_part #extracted products. Originally benzene was the product."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For decomposition of *CGRContainer* back into *ReactionContainer* `CGRpreparer` class can be used. `CGRpreparer` is callable object"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CGRtools import CGRpreparer # import of CGRpreparer\n",
    "preparer = CGRpreparer()         # initialization of CGRpreparer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "decomposed = preparer.decompose(cgr2) # decomposition of CGR2 into reaction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "decomposed # You can see that water absent in products initially was restored. \n",
    "# This is a side-effect of CGR decomposing that could help with reaction balancing. \n",
    "# But balancing using CGR decomposition works correctly only if minor part atoms are lost \n",
    "# but multiplicity and formal charge are saved. In next release electronic state balansing will be added."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r1 # compare with initial reaction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.4 Queries\n",
    "\n",
    "CGRtools supports special objects for Queries. Queries are designed for substructure isomorphism. User can set number of neighbors and hybridization by himself (in molecules they could be calculated but could not be changed).\n",
    "\n",
    "Queries don't have `reset_query_marks` method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CGRtools.containers import*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m10.reset_query_marks()\n",
    "m10 # ether"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "carb = m10.substructure([5,7,8, 2]) # extract projection of carboxyl fragment\n",
    "carb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "q = QueryContainer(carb)  # convert fragment into query\n",
    "print(q) # QueryContainer don't support depiction yet but signatures can be extracted"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CGRs also can be transformed into Query.\n",
    "\n",
    "`QueryCGRContainer` is similar to QueryContainer class for CGRs and has the same API.\n",
    "\n",
    "`QueryCGRContainer` take into account state of atoms and bonds in reactant and product, including neighbors and hybridization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cgr1.reset_query_marks()        # cgr1 is CGRContaier its labels could be calculated \n",
    "cgr_q = QueryCGRContainer(cgr1) # transfrom CGRContainer into QueryCGRContainer\n",
    "print(cgr_q)                    # print out signature of query"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.5. Molecules, CGRs, Reactions construction\n",
    "\n",
    "CGRtools has API for objects construction from scratch.\n",
    "\n",
    "CGR and Molecule has methods `add_atom` and `add_bond` for adding atoms and bonds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CGRtools.containers import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = MoleculeContainer() # new empty molecule\n",
    "\n",
    "m.add_atom('C')  # add Carbon atom using element symbol\n",
    "m.add_atom(6)    # add Carbon atom using element number. {'element': 6} is not valid, but {'element': 'O'} is also acceptable\n",
    "m.add_atom({'element': 'O', 'charge': -1}) # add negatively charged Oxygen atom. Similarly other atomic properties can be set\n",
    "\n",
    "# add_atom has second argument for setting atom number. \n",
    "# If not set, the next integer after the biggest among already created will be used.\n",
    "m.add_atom({'element': 'Na', 'charge': 1}, 4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m.add_bond(1, 2, 1) # add bond with order = 1 between atoms 1 and 2\n",
    "m.add_bond(3, 2, {'order': 1}) # the other possibility to set bond order"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m.calculate2d() #experimental function to calculate atom coordinates. Has number of flaws yet\n",
    "m"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Reactions can be constructed from molecules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r = ReactionContainer() # empty reaction\n",
    "r.reactants.append(m1) # add reactant\n",
    "r.products.append(m11)  # add product\n",
    "# or\n",
    "r = ReactionContainer(reactants=[m1], products=[m11]) # one-step way to construct reaction\n",
    "# or\n",
    "r = ReactionContainer([m1], [m11]) # first list of MoleculeContainers is interpreted as reactants, second one - as products"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# reactants, products, reagents attributes are list-like. \n",
    "r.products.append(m.copy()) # One can add (or remove) molecules directly to this list of products returned by r.products\n",
    "r.flush_cache()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r # coordinates will left unchanged. Thus depiction could look wrong."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r.fix_positions() # this method fixes coordinates of molecules in reaction\n",
    "r"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*QueryContainers* can be constructed in the same way as *MoleculeContainers*.\n",
    "\n",
    "Unlike other containers *QueryContainers* additionally support atoms, neighbors and hybridization lists."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "q = QueryContainer() # creation of empty container\n",
    "q.add_atom({'element': ['N', 'O']}) # add atom list: N or O atom, any isotope and multiplicity, neutral, \n",
    "                                    # number of neighbors and hybridization are irrelevant\n",
    "q.add_atom({'element': 'C', 'neighbors': [2, 3]}) # add carbon atom, any isotope and multiplicity, neutral, \n",
    "                                                  # has 2 or 3 non-hydrogen heighbors\n",
    "q.add_atom({'element': 'O'}) # add oxygen atom, any isotope and multiplicity, neutral, \n",
    "                             # number of neighbors and hybridization are irrelevant \n",
    "q.add_bond(1, 2, 1) # add single bond between atom 1 and 2 \n",
    "q.add_bond(2, 3, 2) # add double bond between atom 1 and 2 \n",
    "# any amide or carboxilic group will fit this query"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(q) # print out signature (SMILES-like)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.6. Extending CGRtools\n",
    "\n",
    "You can easily customize CGRtools for your tasks.  \n",
    "CGRtools is OOP-oriented library with subclassing and inheritance support.\n",
    "\n",
    "As an example, we show how special marks on atoms for ligand donor centers can be added."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CGRtools.containers import MoleculeContainer\n",
    "from CGRtools.attributes import Atom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class MarkedAtom(Atom):  # this class will inherite Atom class\n",
    "    __slots__ = '__mark' # all new attributes should be slotted!\n",
    "    \n",
    "    def __init__(self, **kwargs):\n",
    "        super().__init__(**kwargs)\n",
    "        self.__mark = None # set default value for added attribute\n",
    "\n",
    "    @property\n",
    "    def mark(self):        # created new property \n",
    "        return self.__mark\n",
    "    \n",
    "    @mark.setter\n",
    "    def mark(self, mark):\n",
    "        # do some checks and calculations\n",
    "        self.__mark = mark"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class MarkedMoleculeContainer(MoleculeContainer): # MarkedAtomContainer inherits MoleculeContainer\n",
    "    node_attr_dict_factory = MarkedAtom # override atom container"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = MarkedMoleculeContainer() # create newly developed container MarkedMoleculeContainer\n",
    "m.add_atom('C')               # add atom C using methods inherited from MoleculeContainer\n",
    "m.add_atom('N')               # add atom N\n",
    "m.add_bond(1, 2, 1)           # add bond\n",
    "\n",
    "m.atom(2)['mark'] = 1  # set mark on atom. In this example dictinary setitem supported but update not."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m.atom(2).mark # one can return mark"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "cgrtools",
   "language": "python",
   "name": "cgrtools"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
